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Abstract 

A system of particles is studied in which the stochastic processes are 
one-particle type-change (or one-particle diffusion) and multi-particle an- 
nihilation. It is shown that, if the annihilation rate tends to zero but the 
initial values of the average number of the particles tends to infinity, so 
that the annihilation rate times a certain power of the initial values of 
the average number of the particles remain constant (the double scaling) 
then if the initial state of the system is a multi-Poisson distribution, the 
system always remains in a state of multi-Poisson distribution, but with 
evolving parameters. The large time behavior of the system is also inves- 
tigated. The system exhibits a dynamical phase transition. It is seen that 
for a fc-particle annihilation, if k is larger than a critical value k c , which is 
determined by the type-change rates, then annihilation does not enter the 
relaxation exponent of the system; while for k < k c , it is the annihilation 
(in fact k itself) which determines the relaxation exponent. 
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1 Introduction 



People have studied reaction-diffusion systems, using analytical techniques, ap- 
proximation methods, and simulations. A large fraction of exact results be- 
long to low-dimensional (specially one-dimensional) systems, as solving low- 
dimensional systems should in principle be easier [1-13]. Despite their simplic- 
ity, these systems exhibit rather rich and non-trivial dynamical and stationary 
behavior. Studies on the models far from equilibrium have shown that there 
is a remarkably rich variety of critical phenomena [1]. Among the important 
aspects of reaction-diffusion systems, are the stationary state of the system (or 
one of the quantities describing the system) and the relaxation behavior of the 
system towards this configuration. 

Field theoretic methods have been used to study diffusion, recombination, 
and other dynamic manifestations of non-quantum-mechanical objects (for a 
review see [14]). The techniques of field theory on a lattice are used to examine 
the diffusion and reaction processes of particles [14]. The field theoretic methods 
and the dynamic renormalization group (RG) have been applied to study the 
universal scaling properties of reaction-diffusion models [15,16]. They have also 
been used to study fluctuations in reaction-diffusion problems, for example to 
study the single-species annihilation of k particles to I particles (I < k) [15]. 
In [17] a model is investigated where organisms die and give birth with equal 
rates, and also diffuse. It is shown there that at dimensions smaller than or 
equal to 2, there is aggregation, while at larger dimensions there is no clustering. 
In [18, 19], a model is investigated where particles are created, annihilated, and 
diffused on a lattice. For the case of one- or two-particle annihilation and 
creation, exact results are obtained. In all these cases, bosonic formulations 
have been used, meaning that each site can be occupied by more than one 
particle. 

In this paper, a stochastic model is considered in which the (stochastic) 
variables of the system are the numbers of various types of particles. The term 
type can refer to species as well as position of the particle. There are single- 
particle type-changes, as well as fc-particle annihilations. Specially, a case is 
studied where the annihilation rate tends to zero while this rate times the (k — 
l)'th power of the initial values of the average particle-numbers remain constant. 
It is shown that in this double-scaling limit, the evolution equations for the 
annihilation operators contain only annihilation operators. Specifically, if the 
initial state of the system is a multi-Poisson distribution, then the system always 
remains in a state of multi-Poisson distribution, but with evolving parameters. 
The parameters evolve in a set of equations which in fact are the mean field 
equations. 

It is further shown that the system exhibits a dynamical phase transition. 
The large time behavior of the system is controlled by the spectrum of the 
evolution operator corresponding to single-particle type-changes, and k. It is 
shown that if k exceeds a critical value k c , which is determined by type-change 
rates, annihilation does not enter the relaxation exponent of the system; while 
for k < k c , it is the annihilation which determines the relaxation exponent. 
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The scheme of the paper is the following. In section 2, some general tech- 
niques are introduced, mainly to fix notation. In section 3, the double-scaling is 
discussed. In section 4, the large-time behavior of the system and the dynam- 
ical phase transition are investigated. Section 5 is devoted to the concluding 
remarks. 



2 The general technique 

To fix notation, let's briefly introduce the general technique. Consider a system 
the state of which is characterized by a set of nonnegative integers (say the 
number of different particles). The vector space corresponding to such a system 
can be constructed using the raising and lowering operators (A^'s and A^'s, 
respectively), and the corresponding number operators JV, through 

[A»,A V ] =0, 
[Al At] =0, 
A" |0) =0, 
TV" :=A> Ma ( 3 Al A 13 , 

\n) :=n(4)""|0), 

(m|n) =S mn , (1) 

where the kets \n) form a basis for the vector space corresponding to the system. 
Note that by several types of particles, one may mean several species of particles 
or particles in several places (or both) . 

Corresponding to any set of probabilities P n of finding the system in the 
state \n) (having particles of type fi), there is a probability vector in the 
vector space which is 

\P)=Y,Pn\n). (2) 

n 

Any physical state of such a system is characterized by a vector like 

|P) = /(At)|0), (3) 

where f(A>) is a Taylor series in A' with nonnegative coefficients and with the 
sum of coefficients equal to one. This last condition can be written as one of 
the following equivalent forms 

f(S) =1, 

(S\f(A^)=(S\, (4) 
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where S is a covector all of its coefficients are equal to one, and 

(S| = (0|e s ^. (5) 

The observables of such a system are functions of the number operators. 
The expectation value of the observable g(N) is 

(g(N)) = {S\g(N)\P), 

= (0\g(N + A)f{Al+S)\0), 
=:(0\g{A)f{At+S)\0). (6) 

where in the last equality commutation relations between A's and A^'s have 
been used to rearrange them in g(N + A) so that A''s are a n in the left of 
A's. Specially, if the system has a multi-Poisson probability distribution with 
parameters A M : 

|P)=e A "K-^)| ), (7) 

then 

(g(N))=g(A). (8) 

A general continuous-time stochastic process is described by a linear operator 
(Hamiltonian) H with nonnegative nondiagonal elements and the property that 

(S\ H = 0. (9) 

Such a Hamiltonian can be written in terms of the annihilation and creation 
operators. Specifically, a process involving the annihilation of k particles and 
creation of I particles is described by the Hamiltonian 

H ={Al t ■ ■ ■ Al H A 01 ■ ■ ■ A^ 

- S ai ■ ■ ■ S ai A* •» 4i • • ■ A /3fc i* Sk < ' ' ' 4* A§1 ■ ■ ■ A s «)C^ Pl ... 0k (N), 

=[< • ■ • At, A* ■ ■ ■ A^ - S ai ■ ■ ■ S ai Af(N^ ■ ■ ■ Nf>>)] C a ^ Pl ... Pk (N), 

(10) 

where C"s are nonnegative rates, and M means normal-ordering, that is putting 
At's at the left of A's. 

The evolution of the state vector of the system (\P(t))) is through 

\P(t)) = U(t,0)\P{0)), (11) 

where 

^U(t,0) =HU(t,0), 

f/(0,0)=l. (12) 

So the expectation value of an observable Q at the time t can be written like 

(Q(t))=(S\Q\P(t)), 

={S\Q*(t)\P(0)), (13) 
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where 



Q H (i) :=U- 1 {t,0)QU(t,0), 
^Q H (t)=[Q H (t) > ff H (t)] ) (14) 

One notes that the Heisenberg operators Q H are in fact the ordinary operators 
in them A's and A''s are substituted by A(t)'s and A^(t)'s. 
One also has 

(Q(t)> =(0|Q|P(t)), 

= (0|Q H W|P(0), (15) 



where 



\P) 



e s » A " |P), 



Q H (t) :=U- 1 {t,0)QU(t,0), (16) 



and J7 is defined similar to (|12p. but with _ff in place of H. It is seen that the 
effect of tilde on the operators is just to change A^ a to (A' a + S a ). 

3 Double scaling in annihilation processes of low 
rates 

Consider a reaction-annihilation process with the Hamiltonian 

H = H +H 1 , (17) 



where 



where 



Hq :=M a p A^ a A 13 , (18) 

Hi := J2 G ^-h [A 01 • • • ^ fc - M{N^ ■■■N^)j, (19) 

k 

S a M a a = 0. (20) 



Since only the symmetric part of C enters the Hamiltonian, from now on it is as- 
sumed that C is symmetric. Hq describes a reaction, change of a particle of type 
(3 to a particle of type a with the rate M Q ^, while H\ describes annihilations. 
For the observable g{N), one has then 

(g[N(t)]) = (0\U-\t, 0) g(A) U(t, 0) |P(0)>, (21) 
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where the evolution of U is governed by by H: 

H = M a p Al A? + C Pl ... Pk {A 131 ■■■A^ -AT[(N^ + A^) ■ ■ ■ {N?* + A 13 ")}}, 

k 

(22) 

and 

|P(Q)) = /(At)|0>. (23) 

Suppose that the initial state vector describes a large number of particles, and 
the annihilation rates are small, specifically so that f(A'/X) and (A fc_1 Cp 1 ...p k ) 
(for all fc's) both exist as A — ■> oo. One can then define another pair of annihi- 
lation and creation operators through the transformation 

A 13 =: A a' 3 , 

4=:A-V a . (24) 

Writing A's and A*'s in terms of a's and <r's, and sending A to infinity, it is 
seen that in the Hamiltonian H only those terms survive that are linear in A' . 
So in this double-scaling limit, one can use instead of H the Hamiltonian H s : 

H s := M a pAl t A p - ^kC !3l ... l 3 k N Pl A p °- ■■■A /3k , (25) 

k 

where use has been made of the fact that C"s are symmetric. It is now easily 
seen that in this limit, 

(g[N(t)]) = {0\9[Mt)}\P(0)), (26) 

where A"(<)'s satisfy 

M a p ~^ a pY,kC 01 ...p k A^---A^ if. (27) 

This is a set of differential equations for A" 's, which are commuting at t = and 
remain commuting at later times. Specifically, if the initial state of the system 
is a multi-Poisson distribution (|7|), then with this evolution the system always 
remains in a state of multi-Poisson distribution, but with evolving parameters 
A a (t) which satisfy 

' M a - A^ a p Y^kCp t ...p k A^-.-A^ ) A' 3 . (28) 

This equation can be solved perturbatively. One defines 

A a (t)=:R a p(t)Y^(t), (29) 

where 

d^= M < 

R(0) =1. (30) 
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To be more specific, let's consider a system described by a reaction and a k- 
particlc annihilation. One has 

1 Y a = (R- 1 )^ I>V--/fc (R^ ai Y^)--- (Rf>» ah Y ak ), (31) 

where 

D a fj 1 -p k := -kC v{ p 2 ... h A I/a Plh (32) 

and the (/?2 • • • ftk Pi) means that part which is symmetric with respect to the 
indices. Eq. (|3ip can be rewritten like 

Y"(t) = A«(0)+ f dt' (R-'TAf) D° fil ...p k 
Jo 

x [R* ai (f ) Y^ (f )] ■ ■ ■ [B?> ak (?) Y ak (t% (33) 

and from that 

A a (t) =R a p(t)A< 3 (0)+ f dt'R a a {t-t') D" 01 ...p h K^{t')---K^{t l ). (34) 
Jo 

The above expression can be visualized by a set of graphs. Each graph consists 
of vertices and directed links. Each vertex has one outgoing link and k incoming 
links. Each graph is connected, has no loops, and has only one outgoing link. 
A" (t) is the sum of possible such graphs, the values of them are calculated using 
the following rules. 

• To each point of the graph (the beginning points, the end point, and the 
vertices) is assigned a time. The time corresponding to the end of a link 
should not be smaller than the time corresponding to the beginning of 
that link. The time corresponding to the beginnings of the incoming links 
of the graph are 0, and the time corresponding to the end of the outgoing 
link of the graph is t. 

• To each directed link is assigned a factor R, the argument of which is the 
time corresponding to the end of the link minus the time corresponding 
to the beginning of the link. 

• To each vertex is assigned a factor D. 

• To the beginning point of each incoming link of the graph is assigned a 
factor A(0). 

• the value assigned to a graph is the product of the values assigned to 
various parts of the graph, integrated over the times corresponding to the 
vertices. 

Using this scheme, one can in principle find A(f) up to desired order (number 
of vertices) . 
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4 The large time behavior of the system 



The real parts of the eigenvalues of M are nonpositive, and zero is an eigenvalue 
of M. Assuming that the only eigenvalue with nonnegative real part is zero, and 
that this eigenvalue is nondegenerate, the large time behavior of R is simple. 
The large time behavior of R depends on the spectrum of M near zero. If 
there is a gap in the spectrum, that is if the supermum of the real parts of the 
eigenvalues (apart from zero) is negative, then 

lim R a a{t) =u a Sp, (35) 

t — >oo 

where u is the normalized right eigenvector of M corresponding to the eigenvalue 
0: 



M u =0, 
S a u a =1. 



(36) 



If there is no gap in the spectrum of M near the eigenvalue zero, but still the only 
eigenvalue of M with zero real part is zero, and this eigenvalue is nondegenarte, 
then 

X 

t 



R a (t) 



u a S, 



0' 



t — > OO. 



(37) 



where r and 8 are constants depending on the behavior of the spectrum of M 
near zero. 

In the case there is a gap in the spectrum, for t — ► oo one can substitute the 
right-hand side of ([35]) for R in the expressions for the graphs, as for most of 
the times, the argument of R is large. This is equivalent to rewriting ([34]) as 



A a (t) = u a 



S (3 A?(0)+ / dt'S a DV"AA*(f)-A A (0 



(38) 



A(i)=5^A^(0) 



dt' V [A(0] 



where 
and 



A a (t) =: u a X(t), 



V S a D° Pl ... Pk i/i 
It is easy to solve (l39l) . One has 

dA 



dt 



from which one obtains 
\(t) = 



V \ 



S«A a (0) 



{1 - (jfe - 1) V [Sf3 A /3 (0)] fe_1 i}V(*-i) ' 



(39) 

(40) 
(41) 

(42) 
(43) 
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or 



A a (t) = 



u Q S' 7 AT(0) 



(44) 



{1 + k (k - 1) C 01 ... Pk wft • • • vP* [Sp A/ 3 (0)] fc - 1 t}V(*-i) ■ 

If there is no gap in the spectrum of M near zero, one has to substitute (|37[) 
in (|34p . It is then seen that A(t) is proportional to u, for large times. Putting 
the ansatz 

A Q (t) ~i-"u Q , t^oo, (45) 
along with (1371) in (|34| . one arrives at 

t-S 



t~ 



ct~ 



where 



I(t) 



t-v 



dt' (t-t'y 5 t 



f\ — S ,f—k fi 



(46) 



(47) 



x and y are introduced to ensure that the approximations used for R and A are 
valid in the integration domain. A dimensional analysis shows that for large t, 



So, 



lit) 



t-" 



Cl t- S +C 2 t- kf * 



il—5—k (jl 



est* 



t 



t 



(48) 



^ v -r o 3 v . u^^O. (49) 

The meaning of this, is that the two largest exponents entering this expression 
should be equal. As k > 1, the exponent —k fi is smaller than — fj,. So it cannot 
be among the largest exponents. There remains three possibilities: 



pb = < 5, 



Kk-1 



k5 < 1 



kd > 1 



kS > 1 



(50) 



It can be shown that the third case does not occur. To see this, consider the 
integration corresponding to the a vertex all of its incoming links are incoming 
links of a graph. The time dependence of the integral involved is 



/'(f) - di t- a +d 2 t~ kd +d 3 t 



,1-S-kS 



(51) 



where t is the time corresponding to the end of the outgoing link. It is seen that 
for k5 > 1, the largest exponent in the right-hand side is —S, which shows that 
the result of the integration is proportional to t~ s (for large times). Repeating 
this for successive vertices, One finds that all graphs are proportional to t~ s . 
So the correct value for fi is S. One can then summarize (|50[) in 



fi = max I (5, 



1 



(52) 
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Defining 

k c := J, (53) 

it is seen that the system exhibits a dynamical phase transition: for k > k c , 
annihilation does not enter the relaxation exponent of the system; while for 
k < k c , it is the annihilation which determines the relaxation exponent. 

A note is here in order. If instead of (I45]) , one choose an ansatz that the 
relaxation of A is exponential rather than power law, then the integral on the 
right-hand side of (|34j) tends to zero faster than A itself (as k > 1), which means 
that for large times, only the first term on the right-hand side of (|34|l determines 
A, so it should decay exponentially, which is not the case. 

As an example, consider a system consisting of particles of a single species 
diffusing on a d-dimensional lattice with symmetric rates. Suppose that there is 
a /c-particle annihilation (double-scaled) as well. The types of the particles are 
just the sites of the lattice, denoted by x (d-tuples of integers). The matrix M 
describing the diffusion is then 

d 

M = Y j r l {T l + T- l -2), (54) 
i=i 

where Ti is the one-site translation in the i-th direction. The eigenvalues of M 
are 

d 

E{fi)'.= Y^2n{cas6i-V). (55) 

i=l 

For a finite lattice, 0's are discrete and there is a gap in spectrum at zero. For 
an infinite lattice, the spectrum is continuous at zero, and a steepest-descent 
study shows that 

5= d -. (56) 

This shows that in this case, the system never crosses the critical point k = fc c , 
as k c < 2 and k > 2. 



5 Concluding remarks 

A system was investigated consisting of several types of bosonic particles. By 
bosonic, it is meant there can be more than one particle of each type (at each 
site). In [17-19], similar systems were investigated and exact results including 
phase transitions were obtained for the case of at most two particle creation or 
annihilation. In the case investigated here, there is no creation, but there is 
fc-particlc annihilation, where k can be larger than 2. The case of small anni- 
hilation rate, together with large initial number of particles was investigated in 
more detail. It was shown that this system exhibits a dynamical phase tran- 
sition, which is controlled by k and another parameter related to the rates of 
one-particle reactions. 
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